Avalanches in the lung: A statistical mechanical model 
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We study a statistical mechanical model for the dynamics of lung inflation which incorporates 
recent experimental observations on the opening of individual airways by a cascade or avalanche 
mechanism. Using an exact mapping of the avalanche problem onto percolation on a Cayley tree, 
we analytically derive the exponents describing the size distribution of the first avalanches and test 
the analytical solution by numerical simulations. We find that the tree-like structure of the airways 
together with the simplest assumptions concerning opening threshold pressures of each airway, is 
sufficient to explain the existence of power-law distributions observed experimentally. 
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Recent interactions between physics and physiology 
have resulted in advances in understanding some "sim- 
pler" physiological systems [|l|. In particular, consider- 
able progress has occurred in the general area of statisti- 
cal mechanics and pulmonary physiology [^,^ due most 
likely to the unique tree-like connectivity of the airways 



During a forced exhalation, lungs deflate to very low 
volumes and many peripheral airways close up |5[|. In 
lung disease, closure occurs even during normal breath- 
ing; the closed airways do not reopen for a significant 
portion of the following inhalation |^. As a consequence, 
a large portion of the alveolar space can remain closed 
during the entire breathing cycle leading to severe hy- 
poventilation and an imbalance between ventilation and 
perfusion. The process of opening a single airway is a 
local and isolated phenomenon. However, the dynamics 
of consecutive airway openings in the lung is a highly 
cooperative process. There is recent evidence suggesting 
that during inflation the resistance to airflow of the small 
airways decreases in discrete jumps (7|j^. Thus airways 
do not open individually, but in a sequence of bursts or 
"avalanches" |^ involving many airways; both the sizes 
of these jumps and the time intervals between jumps fol- 
low power law distributions In this paper, we argue 
that the existence of power laws in lung inflation can 
arise directly from the tree-like connectivity of the air- 
ways. We observe that the dynamics of lung inflation 
can be usefully described by a percolation problem on a 
Cayley tree, with the inflated lung volume correspond- 
ing to a percolation cluster. Using this exact mapping, 
we analytically derive the exponents describing the size 
distribution of the first avalanches, and test our results 
using simulations. 

Morphological data [Q show that human (as well 
as other mammalian) lungs constitute an asymmetric 
branching airway structure with « 35 generations. Com- 
plete airway closure appears to occur only in the last 
« 10 — 14 generations where the branching structure 
is largely symmetric [|j . Accordingly, we model this part 



of the airway tree as a binary Cayley tree with airway 
segments that can be either closed or opened (Fig. 1). 
At time i = 0, all airways are assumed to be closed. Lung 
inflation is simulated by applying an external pressure Pe 
at the root of the tree, and gradually increasing Pe at 
a uniform and slow rate. Airways are labeled with 
a generation number i {i — 0,...,N), where N is the 
order of the tree {i — denotes the tree root) , and a col- 
umn number j (j = 0, . . . , 2* — 1). An opening threshold 
pressure Pij is also assigned to each airway (i, j). Exper- 
iments on flexible tube airway models [ pO| conflrm that 
the opening of a single airway is a dynamic process, with 
each airway characterized by a critical pressure threshold 
such that if Pe exceeds this threshold, then the airway 
opens in a short time, which is considered to be instanta- 
neous [ pd] |. Opening of airway {i,j) occurs whenever Pij 
is smaller than the pressure in its parent. 

We assume that Pij is uniformly distributed between 
and 1 ll^, and allow Pe to increase from to 1 in small 
increments. When Pe first exceeds Pqo, the airway (0, 0) 
opens and its pressure is set equal to Pe. Next, the two 
airways (1,0) and (1,1) are tested to see if they can be 
opened with this value of Pe — i.e., if Pe > Pio and/or 
Pe > Pii- If one or both conditions are met, then the 
airways (1,0) and/or (1, 1) are also opened. This open- 
ing is then continued sequentially down the tree until no 
airway is found with its Pij < Pe. 

Of particular interest is the fact that a small increase 
in Pe can lead to an "avalanche" in which many airways 
open simultaneously. When the first avalanche stops, Pe 
is further incremented and pressures in the open airways 
are updated. We iterate this process until all airways 
open. The location and size of the next avalanche de- 
pends on the distribution of Pij in the accessible region. 

We do not treat the full problem analytically, but we 
can obtain exact results for the distribution of the first 
avalanche. At f = 0, we increase Pe until the first 
avalanche occurs and we calculate its size s. Then we 
restart the simulation with a new set of thresholds {Pij}- 

Before we consider two possible definitions of s, we note 
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that gas exchange in the lung occurs only in the "opened" 
alveoli (the terminal units of the bronchial tree) which are 
in communication with the trachea. For this reason, in 
Definition A, s denotes the number of alveoli, defined as 
the number of elements in the last generation, TV, that 
become open. Motivated by percolation theory |]l3| , |l4| , 
in Definition B, s is the number airways that open fol- 
lowing an increase of Pe that opens at least one airway. 
The physiological rationale for definition B is that when 
the lung is deflated to low volumes, most airways close. 
However, often there remains trapped air in the alveoli. 
Thus concerning gas exchange it may not be necessary 
that an avalanche reach the bottom of the tree for it to 
connect alveoli with the trachea. 

We study n(s), the size distribution of first avalanches. 
For Definition A, 11'^ (s) shows a single power law behav- 
ior with an exponent ■y^ = 0.9 {~ I.O) (Fig. 2). For 
Definition B the function n^(s) has two regions (Fig. 2): 
a first region with a steep power-law decay and a second 
region with a moderate power-law decay, with a crossover 
at a scale TV, 



n^(.) 



[s < N] 
[s>iV] 



(1) 



The exponent ~ 1.9 (« 2.0) for the first regime, while 

7^ = 0.9 (« 1.0) for the second regime which extends 

to sizes including all branches, i.e., almost to a size of 
2N+1 _ ^ 

We argue that for Definition B, this problem can be 
mapped onto the percolation problem for the Cayley tree 
I p^ . In the percolation problem, we occupy randomly 
every branch of the tree with a probability p. Then, 
starting from the root, we connect all occupied branches 
that are neighbors of each other. Definition B concerns 
the cluster of connected bonds that starts from the root. 
The size of this cluster depends on the fraction p of oc- 
cupied bonds. As we approach a critical probability Pc, 
the typical size of a cluster can be characterized by the 
Styp ~ \p — pc\~^^'^- Both a and pc can be calculated ex- 
actly due to the branching nature of the tree: — \ and 
= i [Q. In general, the size distribution of the finite 
clusters in the infinite system obeys the scaling form fl^ ] 



n{s)^s--f{s-\p-p,\), 



(2) 



where t ~ 3/2, and /(u) = const for u ^ 1 and 
/(w 1) — » 0. To connect percolation theory to the 
lung model, instead of occupying randomly the branches 
with probability p, we assign a random number or pres- 
sure threshold value to each airway. We then define a 
cluster to be the set of airways that have a threshold 
smaller than a predefined value p and are connected to 
the root. When Pe exceeds Pqo, we open all airways be- 
low the root which have a threshold value smaller than 
-Poo- 

If Pqo is fixed and set equal to p, then this is exactly 
the percolation problem on the Cayley tree and the dis- 
tribution of the cluster sizes or avalanches is given by 



(g). However, in our case Pqo is also a random variable. 
Thus, in order to obtain the size distribution of the first 
avalanche, we must integrate the cluster distribution over 
the probability p, from to 1 with the result that 
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which predicts 7^ = 2, in agreement with the scaling 
observed for s < in Fig. 2. 

These calculations assume that the system size is in- 
finite. No avalanche with s < N can reach the bottom 
of the tree, so the scaling behavior for s < A is that of 
the infinite tree with an exponent = 2. On the other 
hand, avalanches of size s > N are affected by finite size 
effects, and indeed the data for s ^ N indicate a different 
exponent. Moreover, finite-size effects will always affect 
the scaling behavior of H'^(s), since every avalanche that 
leads to the opening of one or more alveoli must open at 
least A^ airways. 

Percolation theory on the Cayley tree does not pro- 
vide the scaling behavior of the cluster size distribution 
dominated by finite size effects. Nevertheless, we can 
obtain the scaling exponent using a generating function 
approach, developed in the theory of branching processes 
. We define the generating function of order A" 
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where P^'^ {p, s) give, for definitions A and B, the prob- 
ability that in a tree with A^ generations we have an 
avalanche of size s for a given Pqo = p. Therefore, 



Hat (s) 



f P^ ''{p,s)dp. (5) 
Jo 

These generating functions satisfy the recursion relations 

.g^+ifea;) = x'''[{l - p) + pg^'^ {p,x)f , g^'^{p,x) = x 

(6) 

where cj) — Q for definition A and cj) — 1 for Definition 
B, g^ is a polynomial in x of degree 2^ and gfj is a 
polynomial of degree 2^+^ — I. We obtain the distribu- 
tion functions H'^(s) and H^(s) by numerical integration 
of the coefficients of these polynomials with respect to 
P — Poo- The results for H"*(s) and H^(s) are shown 
in Fig. 2; note the good agreement between simulations 
and theory, despite the fact that there are no adjustable 
parameters in the calculation (the theoretical line being 
determined solely by the value of A). 

For A — i- 00, for any x < 1, the generating function 
gfj{p,x) approaches the limit g^{p,x) = [l — 2p{l—p)x — 
\/l ~ 4p(l — p)x ]/2p'^x, which can be expanded in pow- 
ers of X. On integrating the coefficients of this expan- 
sion with respect to p we obtain H^(s) = l/s(2s -I- 1) for 
s < N, which implies an asymptotic exponent 7;P = 2. 
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Next we consider Definition A, and show that Il^{s) ~ 
1/s for large s, so that 72^ = 1- For large N and large s, 
it follows from general theorems that for p > 1/ over2 



-C{s/s^)M, so^i2pf. (7) 



Here Sq is the average avalanche size (number of open 
alveoli) in the generation N, j{p) is a continuous function 
of p for p < 1, and C is a positive number with a weak 
dependence on (s/sq)- For p < ^, so decays exponen- 
tially with N. Thus, for p < \ the coefficients P^{p, s) 
for large s become negligibly small and do not contribute 
to the integral (H). In contrast, for p > ^ the probability 
of a non-zero avalanche in definition A — which is equal 
to the sum of all the P^{p, s) with s > 1 — is finite when 
N ^ 00 and equal to {2p — This quantity should 

be used as normalizing constant in the equation (|7|). In- 
tegrating equation (|^) with respect to p from ^ to 1 with 
the help of the saddle point approximation, we find 
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n^(s) ^ (1 _ s-i/JV) [Ar«s«2^], 



sN 



(8) 



so 7^ 



1. If we expand Eq. for small s, we find 
(s) ~ lns/(siV^); hence we expect to find an effec- 
tive exponent that is smaller than the asymptotic value 
7^ = 1, and indeed our simulations give 7^ — 0.9 (Fig. 2). 
For very large s, comparable with 2^, the saddle point 
approximation is no longer valid and we observe (Fig. 2) 
the "kink" near the end of the distribution Eq. (i) 
also holds for definition A, so 7^ = 1. Note that Eq. (S) 
is valid for trees with any coordination number. 

Having derived the above exponents analytically, we 
next examine their "universality" by discussing how de- 
viations from the assumptions made in the model may 
affect the scaling exponents, (i) The first assumption 
(which matters only for Definition B) is that we neglect 
the fact that in the lung the length £ and the radius 
r and hence the volume of the airways depend on the 
generation number i Previously, we modeled this 
generational dependence such that £i+i ~ £i/0-9, and 
r^+i — ri/0.86 where the scaling factors (0.9 and 
0.86) arise from morphological data This exponen- 
tial dependence should not affect the scaling behavior, 
an expectation we verified by simulations, (ii) The sec- 
ond assumption, that the distribution of Pij is uniform, 
matters for both definitions. Unfortunately, direct ex- 
perimental data on the distribution of Pij in the lung 
are not available. However, even if the distribution is 
not uniform, but normal or exponential, the scaling ex- 
ponents will not be influenced as long as the values of 
Pij are not correlated. Correlations among Pij have not 
been reported, (iii) While the assumption of a uniform 
distribution is physiologically reasonable, it is also possi- 
ble that there is a weak generational dependence of Pij 
| p2t which can reduce the scaling region and/or change 
the value of the exponents. A stronger generational de- 
pendence of Pij in which the mean of Pij as a function 



of i increases from the root to the bottom by at least a 
factor of 10 will, however, break down the scaling behav- 
ior ijl^ . As a consequence, the very existence of scaling 
exponents found in experimental data provides indi- 
rect evidence that the distribution of Pij does not have 
any generational dependence. 

In summary, we have studied a statistical mechanical 
model of the distribution of the first avalanches during 
lung inflation. Our main result is an analytically-soluble 
model which, compared to the more realistic model of 
Ref. Q, permits exact calculation of the scaling expo- 
nents with the avalanche size defined either as the num- 
ber of alveoli (Definition A) . or number of opened airways 
(Definition B) . We have found that the tree- like structure 
of the airways with the simplest assumptions concerning 
opening theshhold pressures is sufficient to explain the 
existence of power-law distributions observed experimen- 
tally Q . Finally, the fact that the size distribution of the 
first avalanches follows a power law suggests that high 
pressures for at least short periods of inspiration might 
be necessary to open up larger alveolar volumes. Thus, 
our results may also find important applications in the 
design of appropriate waveforms for artificial ventilation 
of patients who suffer from substantial airway closure and 
alveolar collapse. 
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FIG. 1. Schematic diagram of the airways represented by 
a branching tree. The airways are labeled by a generation 
number (i — 0, A*') and a column number {j = 0, 2' — 1). 
A threshold pressure < Pij < 1 chosen from a uniform 
distribution is also assigned to each airway. 

FIG. 2. Double logarithmic plot of the avalanche size dis- 
tributions n(s), obtained by computer simulation on a Cay- 
ley tree of 12 generations. Shown are data obtained for 10* 
realizations for each of the two avalanche size definitions dis- 
cussed in the text. Definition A (closed circles) and Definition 
B (open circles). Also shown, for comparison, are the exact 
results obtained using the generating function approach de- 
scribed in the text. Again, both definitions are shown: Defi- 
nition A (solid line) and Definition B (dotted line). 
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